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4  . 


SUMMARY 


A  computational  method  for  unsteady  aerodynamics  of  fighter 
configurations  at  high  angles  of  attack  is  developed.  The  leading- 
edge  vortices  are  represented  by  free  vortex  filaments  which  are 
adjusted  iteratively  to  satisfy  the  force-free  condition.  The 
small-disturbance,  unsteady  potential  equation  is  solved  in 
frequency  domain  for  motions  in  pitching,  plunging,  flapping,  side 
movement,  rolling,  and  yawing  oscillation  in  compressible  flow. 
Computed  results  in  rolling  moment  coefficients  due  to  side 
acceleration  are  compared  with  data  for  60-deg  and  80-deg  delta 
wings.  Lateral-directional  characteristics  for  an  F-106  B 
configuration  are  also  compared  with  data  obtained  in  forced- 
oscillation  test.  It  is  shown  that  reasonable  results  can  be 
obtained  by  the  present  unsteady-flow  method,  but  not  by  the  steady- 
flow  theory.  Calculation  of  dynamic  stall  effect  on  a  rectangular 
wing  of  aspect  ratio  4  is  demonstrated  by  using  experimental 
sectional  data.  Although  no  data  for  the  wing  are  available  for 
comparison,  the  results  appear  to  be  plausible. 
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Dimension 


Freestream 


1 .  INTRODUCTION 


When  a  fighter  aircraft  performs  maneuvers  at  high  angles  of 
attack,  the  aerodynamic  forces  and  moments  acting  on  the  aircraft 
are,  in  general,  dependent  nonlinearly,  not  only  on  the  motion 
variables  but  also  on  the  time  rate  of  change  of  these  motion 
variables.  To  predict  the  dependency  on  motion  variables,  such  as 
the  dihedral  effect  (C  ),  the  VORSTAB-II  code  was  developed  (ref. 

3 

1).  The  code  is  based  on  a  steady-flow  theory,  so  that  the  unsteady 
aerodynamic  effect  on  stability  parameters  cannot  be  accounted 
for.  This  shortcoming  is  particularly  serious  in  predicting  dynamic 
stability  parameters,  such  as  roll-damping  derivatives,  at  high 
angles  of  attack  (ref.  1).  This  is  because  in  an  unsteady  rolling 
motion  about  the  body  axis,  dynamic  change  in  sideslip  will  be 
induced  to  result  in  rolling  moment  coefficient  due  to  side 
acceleration  (i.e.  C  ).  The  latter  may  provide  additional  damping 

3 

or  undamping,  depending  on  configurations.  In  fact,  significant 
effect  of  C  on  flight  dynamics  has  already  been  illustrated  in  the 

3 

past  (refs.  2  and  3).  These  dynamic  stability  parameters  can  only 
be  calculated  through  the  unsteady  aerodynamic  theory. 

In  general,  the  unsteady  aerodynamic  theory,  relative  to  a 
sttady-flow  theory,  provides  the  following  mechanisms  to  affect  the 
dynamic  stability  parameters: 

•  aerodynamic  lag  in  attached  flow; 

•  vortex  lag  in  vortex  flow,  including  dynamic  vortex- 
breakdown  effect; 

1 


•  dynamic  boundary  layer  and  stall  effect  in  viscous  flow. 
The  aerodynamic  lag  in  attached  flow  is  produced  by  convective 
effect  of  shed  vortices  and  the  finite  speed  of  disturbance 
propagation  in  a  compressible  flow.  This  effect  is  automatically 
included  in  any  appropriate  inviscid  unsteady  aerodynamic  theory  in 
compressible  flow. 

On  the  other  hand,  vortex  lag  effect  is  important  when  a  wing 
with  edge-separated  vortex  flow  is  in  an  unsteady  motion.  This 
effect  has  been  modeled  in  the  past  through  the  method  of  suction 
analogy  in  compressible  flow  (ref.  4),  the  method  of  separated 
vortex  filaments  in  incompressible  flow  (refs.  5-7),  and  the 
solution  of  Navier-Stokes  equations.  The  method  of  references  5-7 
has  been  illustrated  only  for  simple  wing  planforms  in 
incompressible  flow  and  is  based  on  a  time-step  integration  for  the 
solution.  The  solution  of  Navier-Stokes  equations  for  airplane 
configurations  is  too  "computation  intensive"  for  applications  to 
preliminary  design. 

Several  empirical  modeling  methods  for  dynamic  stall  effect  on 
airfoil  sections  have  been  proposed  as  summarized  in  reference  8. 
Accuracy  of  Navier-Stokes  solutions  in  two-  and  three-dimensional 
flows  is  uncertain,  in  particular  in  deep  dynamic  stall  regime. 

In  this  investigation,  the  effect  of  vortex  lag  is  modeled  by 
using  unsteady  discrete  vortex  filaments  in  compressible  flow.  This 
is  done  by  converting  the  VORSTAB-II  code  into  an  unsteady  version 
in  frequency  domain.  The  main  objective  is  to  calculate  dynamic 


stability  parameters.  For  the  dynamic  stall  effect  on  a  wing,  a 
method  based  on  using  2-D  test  data  is  proposed. 
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2.  THEORETICAL  APPROACH 

2.1  General  Concept 

The  basic  methodology  considered  in  the  present  investigation 
is  based  on  the  subsonic  unsteady  potential  equation  with 
correct'  .ns  for  steady  boundary  layer  separation  and  embedded 
oscillatory  free  vortex  filaments  for  vortex  flow  effects.  The 
unsteady  potential  flow  calculation  is  based  on  the  small 
disturbance  equation.  Only  harmonic  motions  will  be  considered.  By 
the  assumption  of  small  amplitudes  of  motion,  the  solution  can  be 
expressed  as  the  sum  of  two  components,  one  representing  steady  flow 
and  the  other  being  oscillatory.  It  is  assumed  that  the  solution 
for  some  sinusoidal  motion  of  a  wing  with  frequency  to  can  be 
expressed  in  the  form  of 

<)>(x,  y,  z,  t)  =  <j>s (x,  y,  2)  +  ^(x,  y,  z)e1U>t  (1) 

where  <J>  is  the  velocity  potential,  and  ij>g  is  the  steady-state 
solution,  which  is  a  function  of  the  spatial  coordinates  x,  y,  and 
z.  ^  is  an  unsteady  component  and  is  a  complex  function,  thus 
allowing  for  phase  differences  between  the  motion  and  the  resulting 
flow  field. 

In  the  present  method  the  lifting  surface  is  represented  by  a 
distribution  of  oscillatory  horseshoe  vortex,  which  is  the  exact 
solution  of  the  unsteady  small  disturbance  equation.  The  strength 
of  the  vortex  can  be  obtained  by  satisfying  the  boundary  condition 
on  the  lifting  surfaces.  To  properly  account  for  the  leading-edge 
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singularity  of  pressure  loading  in  the  linear  theory,  and  hence  the 
leading-edge  thrust,  the  unsteady  Quasi-Vortex-Lattice  Method 
(UQVLM)  is  used  (ref.  4). 

In  addition,  the  oscillatory  vortex  flow  arising  from 
separation  along  the  leading  and  side  edges  of  a  low  aspect  ratio 
wing  is  simulated  by  using  discrete  free  oscillatory  vortex 
filaments.  The  mathematical  expression  for  the  latter  will  be 
derived  later. 

Note  that  unsteady  fuselage  effect  is  not  investigated  in  this 
study. 


2.2  Unsteady  Potential  Flow  Theory 


The  unsteady  small  disturbance  potential  equation  in  a  moving 
coordinate  system  is 


2  2  2  2  2  2 
U4  a  <(>  a  <J>  d  b  1  d  0  2U  d4(|> 

2;  .  2  .2  .  2  2  2  "  2  dx5t 

a  5x  oy  02  a  ot  a 

cts  *  no  co 


0 


(2) 


By  the  application  of  the  method  found  in  reference  4,  the  above 
equation  can  be  reduced  into  the  following  integral  equation: 

T+X 


0(x ,y ,2 , t )  =  -^  //  tfc  / 


r^Tj-  Ap(£  ,T1  ,t - rr-^-)(^-)dT]OSar|) 

MR-x  ® 

0 


1  -  M' 


(3) 


where 


,  2  ^  2  ^  2 
r  =*  /t  +  y  +z 
o  o 


0  is  actually  the  solution  of  equation  (2)  and  represents  a  doublet 
or  oscillatory  horseshoe  vortex  with  the  strength  of  AC^.  By 


5 


introducing  the  harmonic  time  variation  such  that  for  any  t 

iwt 

Ap(5,  n»  tQ)  =  Re[Ap(?,  n)e  °] 

icjt 

4>(S»  h>  tQ)  =  Re[$(x,  y,  z )e  °] 

After  substitution  of  (4)  into  (3): 

0(x,  Y,  z)  =  ff  ^  <(xq>  y,  z,  K,M)d^dTi 


(4) 


(5) 


where  <  is  defined  as 


0  r 

K  dc  J 


+MR 


-  dx 
r 


1-M 

Since  equation  (3)  is  equivalent  to  a  doublet,  the  solution  of 
equation  (2)  can  be  represented  by  doublet  distributions.  The  thin 
wing  approximation  is  used  throughout. 

To  calculate  the  induced  velocity,  due  to  a  unit  strength  of 
doublet  at  any  (x,  y,  z)  location,  the  method  of  reference  4  will  be 
followed.  The  resulting  equations  are 

'  /  C  AC  ( 

J  P 


5<p 

dz 


5<j>  = 
dy 

d<}> 

dx  L  1  p 


x 

f ' 

‘x 

X, 


3<t>1 

a<t>2 

'dz 

dz  ' 

dd1 

d<D2i 

dy 

+  dy  ‘ 

_dd  i 

d<J>2 

'dx 

+  dx 

(6) 


where 
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1  j'  ** 

'l  *8?  l  [<T  +  ^2)Zoe:i!,I‘i"(ulr1  +  x0)/v.)]dn 

L  r,  Rrt 


i(^>  y»  2)  —  -  v  J  2  y,  T| )  dT) 

®  L 


I(x,  y,  z,  ti)  =  /  (1 - 5-  -  ,  )exp  [-iu)(t’  +  x  )/V  ]dt’  , 

Vi  <^'2-rJ)1/2  »  - 

I  denotes  the  summation  overall  spanwise  strips,  and  x^  and  xt  are 
the  x  coordinates  of  the  leading  and  trailing  edges,  respectively, 
of  the  chord  through  the  collocation  (or  control)  points. 

With  the  transformation 

£  =x£  +c(l  -  cos9)/2 

equation  (6)  becomes 


2  J  ' 

o  u 


3d)  c  *  S  (J)  .  S  0  « 

af  (x>  y*  z)  “It/  ^cn  (e)sine(^i-  +  ^)de 


Note  that  sin(9)  cancels  the  square  root  singularities  of  AC  at 

Pu 

the  leading  and  trailing  edges.  Therefore  equation  (5)  can  be 
reduced  to  a  finite  sum  through  the  midpoint  trapezoidal  rule. 

Hence , 

|±  (x.,..)  -  If-  J  AC  0 

etc.  c  k-1  u 

where  Nc  is  the  number  of  integration  points  and 

9  =  (2k  -  1)*/(2N  ). 

x  c 

The  location  of  the  "bounded"  element  of  the  oscillating  horseshoe 
vortices  is  given  by 


^k  =  xi  +  T  (  1  ”  cos  [  ( 2k  -  l)u/2Nj},  k  *  1,  ...,  N_ 
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The  x  coordinates  of  end  points  of  the  bound  elements  are  given 


by 


Xlk  ^lk’  X2k  ?2k 


(9) 


To  eliminate  the  Cauchy  singularity  in  the  chordwise  integral 

(ref.  4),  a  special  set  of  control  points  is  given: 

x  =  x0  +  c[l  -  cos(iu/N  )]/2  i  =  1,  N 

c  x.  c  c 

yj  =  b/2{ 1  -  cos[Jx/(Ns  +  l)]}/2  j  =  1,  N.  (10) 

where  Ng  is  the  number  of  spanwise  strips. 


Boundary  Conditions 


The  boundary  conditions  for  equation  (2)  are  as  follows:  at  an 
infinite  distance  from  the  surface  disturbance  potential  should  go 
to  zero;  and  on  the  wing  surface,  the  normal  velocity  component  to 
the  wing  surface  should  be  zero.  Assuming  that  the  wing  surface  can 
be  described  as 


Z  =  z(x,  y,  t) 

where  Z  is  height  above  the  plane,  (Z  =0),  of  the  wing  surface. 
Assume  that  during  oscillation  every  point  of  the  wing  moves  in  the 
vertical  direction  only.  Thus,  the  velocity  vector  of  a  point  on 
the  wing  is 
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The  velocity  vector  of  the  fluid  on  the  wing  surface  is 

w0=(U+u)i+vj+wk  (12) 

If  these  velocity  vectors  are  resolved  into  tangential  and 
normal  components  on  the  wing  surface,  and  by  applying  the  flow 
tangency  condition,  the  following  equation  will  be  obtained: 


n  = 


dz  ?  dz  2  T 

-  -r—  1  ~  -r—  2  +  k 
dx _ Sy  _ 

A  *  <i>2  +  W 


-*•  *  ■* 
w^  •  n  =  w^  •  n 


(13) 


(14) 


by  substituting  (11),  (12),  and  (13)  into  (14), 

dz  tI  dz  5z 

W  =  -r—  +  U  "5 —  +  V  t— 

ot  <=  ox  oy 

a.  Longitudinal  Motion 

For  a  wing  in  plunging  motion  with  displacement  h(y,  t),  in 
pitching  with  angular  displacement  9(y,  t)  about  x  ®  xa,  and  in 
rolling  oscillations  with  roll  angle  0r(t),  the  total  vertical 
displacement  Z(x,  y,  t)  is  given  by 

Z(x,  y,  t)  =  -h(y,  t)  -  9(y,  t)(x  -  x  )  -  y<j>  (t)  (15) 

a  r 

dz 

Note  that  v  —  can  be  ignored  because  of  small  perturbation 
assumption. 


=  -h  -  9  (x  -  x  )  -  y0 
ot  a 


5z 

3x 


-9(y,  t) 

Hence,  the  nondimensional  normal  velocity  is  given  by 
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h  0  ,  \  7  Z 

W  =  -  —  -  y-  (x  -  xa)  -  70  -  0 


(16) 


Assuming  the  harmonic  time  variation  such  that 

w(x,  y,  t)  =  Re [w(x,  y)eiujt],  etc.  , 
equation  (16)  becomes 

-w(x,  y)  “  Aa  =  p  h  +  p  8(x  -  xfl)  +  i  yt  +  9  (1?) 

r  r  r 

The  additional  angle  of  attack  (Aa)  is  added  to  the  steady 
angle  of  attack  to  produce  the  total  normal  velocity,  (w^_)  ,  which 
is  (ref.  4) 

azc 

w  (x,  y)  =  -r —  cosa  -  sin(a  +  Aa)  (18) 

t  *  J  dx  s  s 

dz 

=  — —  cosa  -  sina  -  Aacosa 
ox  s  s  s 

b.  Lateral-Directional  Motion: 

For  a  wing  in  sideway  motion  with  displacement  d(t),  and  in 
yawing  with  angular  displacement  <f>(0  t)  about  q  =  x  ,  the  total 
side  displacement  Y  is  given  by 

F(x,  y,  z,  t)  =  Y  -  f(x,  z,  t),  (19) 

where 

Y  =  -d(t)  -  4>(t)  (x  -  xa) 


Noting  that 


9f  x  ri  af  - 
dt  ®  dx  V’ 


(20) 


the  following  side  velocity  contribution  from  the  motion  of  aircraft 
will  be  obtained: 


09  09  CO 


(21) 
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If  the  wing  has  dihedral,  the  following  normal  velocity  will  appear 
in  boundary  conditions;  otherwise,  it  will  be  zero. 


w  v  d  a' 

—  cost  -  —  sinT  =  -ik[j—  sinT  + - j -  sinT]  -  <0sinT  (22) 

CD  GO  ■£  2T 

where  t  is  dihedral  angle. 

The  total  normal  velocity  can  be  canceled  on  the  wing  by  using 
a  steady  horseshoe  vortex  distribution  for  the  first  two  terms  on 
the  right-hand  side  of  equation  (18),  and  by  a  distribution  of 
doublet  for  the  last  term  in  equations  (18)  and  (22). 


2.3  Free  Vortex  Filaments 


A  method  similar  to  steady  free  vortex  filament  is  used  to 
model  the  unsteady  leading-edge  vortex  separation.  Here  each  vortex 
filament  will  be  divided  into  many  segments,  and  each  segment  will 
be  called  an  oscillatory  line  vortex.  The  induced  velocity  due  to 
an  oscillatory  line  vortex  of  unit  strength  is  derived  in  Appendix 
A.  The  resulting  induced  velocity  components  due  to  an  oscillatory 
line  vortex  are  given  as  follows: 


5  <{> 

50L 

50, 

+  - 

5x 

5x 

2)x 

50 

50l 

50 

5z 

5z 

5z 

50 

II 

o> 

-e* 

►— 

50 

5y 

5y 

+  dy 

where : 


(23) 
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Ac  -i— (x  r.-Hc  ) 

♦l  =8/f(Fl  +  F3  -  Ve 


iu>  ,  ,  . 

.  -  -rr-  (T,  r  -he  ) 

(~—)F  ^  (Vl  +  Xo)e  dTl1 

L  ®  ' 


f  =  f(+j  +  |)  2  dn 

ri  Rri 


Z  X 

Fi-/4rd” 

ri 


=  Arctan  [- 


z(x1  -  x) 


(yL  “  y)/(x1~x)2  +  02(y1  -  y)2  +  p2z2 


F^  =  Arctan  [- 


z(x2  -  x) 


(y2  -  y W(x2  ~  x^2  +  P2(y2  "  f)2  +  32z2 


>2  *  *  ir  /  S  I(*’  y'  5'  I’)d'1 

L  ri 


where : 


I  =  /  (1 - 5 - 5-775-) exp  [-iu)(t  +  x  )/V  ]dx' 

-r  xr  L  (t'  +  rL) 

As  mentioned  before,  sraall-amplitude  harmonic  motion  is 
assumed.  Therefore,  the  effects  of  oscillatory  free  vortex 
filaments  and  steady  state  free  vortex  filaments  can  be 
superimposed.  In  both  cases  the  force-free  condition  on  vortex 
filaments  will  be  enforced.  Also,  partial  separation  along  the 
leading  edge  is  allowed  (ref.  1). 
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Note  that  the  derivation  of  equation  (23)  is  based  on  the 
assumption  that  the  vortex  line  is  situated  on  the  plane  of  z  =  0. 
Therefore,  the  induced  velocity  is  always  calculated  first  with 
respect  to  a  local  plane  of  reference  containing  the  vortex  line; 
and  the  results  are  then  transferred  to  the  fixed  plane  of 
reference. 

2.4  Method  of  Solution  with  Vortex  Flow 

The  problem  is  resolved  into  steady  and  unsteady  components, 
through  harmonic  assumption.  Then  after  a  known  steady  state 
solution  is  obtained,  the  solution  for  the  unsteady  motion  will  be 
calculated  and  superimposed  to  the  mean  steady  flow.  The  basic 
unknowns  of  the  problem  are  the  strength  of  oscillatory  horseshoe 
vortices,  strength  of  oscillatory  leading-edge  vortices,  and  the 
locations  of  the  free  elements.  The  problem  is  nonlinear  because 
the  locations  of  the  leading-edge  vortices  are  unknown  a  priori. 
Therefore,  the  solution  procedures  based  on  an  iterative  method 
described  below  must  be  used: 

a.  Prescribe  the  doublet  lattice  for  the  wing  surface,  and 
the  steady-state  locations  of  the  free  elements  over  the 
wing. 

b.  By  satisfying  flow  tangency  condition  on  the  wing  surface, 
obtain  the  doublet  and  oscillatory  leading-edge  vortex 
strengths. 

c.  Calculate  all  the  aerodynamic  characteristics. 
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d.  Enforce  force-free  condition  on  free  elements  by  aligning 
them  in  the  direction  of  local  velocity  vector. 

e.  Repeat  steps  b  through  d  until  a  converged  solution  is 
obtained. 

Note  that  typically  the  out-of-phase  induced  velocity  on  free  vortex 

filaments  is  quite  small  compared  to  the  free  stream. 

Therefore,  C  is  calculated  with  the  symmetric  vortex  shape.  In 

0 

addition,  in  the  vortex-breakdown  region,  it  is  assumed  that  the 
phase  angle  is  90  deg,  so  that  the  in-phase  force  component  without 
breakdown  effect  becomes  out  of  phase,  and  vice  versa. 


2.5  Calculation  of  Aerodynamic  Characteristics 


The  pressure  distribution  on  lifting  surfaces  Is  calculated  by 
applying  unsteady  compressible  Bernoulli's  equation  (ref.  9),  which 
has  the  following  form: 


5$  ,  7<0  •  V$  ,  rp  dpl  U°° 

at  2  J  Pl  2 

Po,  1 


(24) 


By  using  the  isentropic  relation,  definition  of  pressure 
coefficient,  and  finally  linearizing  the  equation,  the  following 
relationship  between  pressure  coefficient  and  velocity  potential  is 
derived  (ref.  9): 

2  5$  2  5$  /,c 

Cp  *  ~  U  5x  2  at  ( 

oa  U 

00 

To  evaluate  change  in  pressure  coefficient,  the  relationship 
between  velocity  potential  and  circulation  is  used;  then  the 


14 


resulting  equation  in  terms  of  vortex  strength  T  becomes  (ref.  9,  p. 
210) 


2  rarA  IT  ar, 

Ac  =  +  — o  LT7  + 

p  ^2  at  ®  ox 


(26) 


To  solve  equation  (26),  harmonic  approach  is  used,  to  transform  the 
above  partial  differential  equation  (26)  into  a  first-order  ordinary 
nonhomogenous  differential  equation,  in  which  the  circulation 
amplitude  is  an  unknown,  i.e.: 


2  ftr 

Ac  =  —r  [iuT  +  U  —  i 
p  ^2  ®  5x 


(27) 


The  solution  of  Equation  (27)  is 


r  =  —  e 


1<JX  _  16JX 

U_  U  U 

CD  A  CO 

J  Ac  e  dx 

P 


(28) 


Finally,  the  lifting  pressure  due  to  oscillating  bound  elements 
and  chordwise  vortex  elements,  based  upon  the  following 
relationship,  will  be  obtained. 


(Ac  )B  -  2(u  -  v  tan*)  ~ 

r  CO  CD 


/.  ,  _  lr  v  2ik  T 

(Ac  )T  ^  u  a  u 

CO  J*  CO 


(29) 


(30) 


Note  that  the  calculated  Ac^  is  interpolated  if  necessary  to  obtain 
ACp  at  integration  stations.  The  calculated  pressure  forces  are 
then  used  to  determine  the  aerodynamic  characteristics  in  a  way 
similar  to  the  steadv  flow. 
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2.6  Longitudinal  Dynamic  Stability  Derivatives 


In  general,  the  amplitudes  of  lift  and  pitching  moment 
coefficients  in  oscillatory  motion  can  be  expressed  in  terms  of 
stability  derivatives  as  follows  (ref.  10) 

S  ' ik  r  S,  -  k\  f +  +  ikec», +  •  •  • 

r  h  hr  9  9 

C  =  ik  -  k2C  +  9C  +  ik9C  +  .  .  . 

m  i  C  ml  jL  mn  ms 

r  m^  nr  9  9 

At  zero  frequency,  C  =  C  ;  and  at  small  frequencies, 

m.  m 
9  a. 

C  =  C  +  C 
m*  m  m* 

9  q  a 


where  C  can  be  shown  to  be  equal  to  C  .. . 

% 


(31) 

(32) 


(33) 


2.7  Lateral  Dynandc  Stability  Derivatives 


In  a  similar  approach  as  in  the  longitudinal  case,  the 
amplitudes  of  side  forces,  and  rolling  and  yawing  moment 
coefficients,  in  each  mode  of  oscillatory  motion  are  expressed  in 
terms  of  stability  derivatives.  The  results  are  as  follows: 


a.  Sideslipping  (g  =  d) 


C  =  ik(C  -  k2  )d  -  k2 

y  ya  c  ..  c 
3  ye 


Cn  =  ik(C  ..  )d  -  k  dC 

"(3  n3 


C  =  ik(C  -  k2C  )d  -  k2dc 


(34) 


(35) 

(36) 
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b.  Rolling  Motion  (p  =  $) 


C  —  ikij)  (C 

-  k2C  ..) 

y  y 

y„ 

P 

p 

C  =  ik4>(C 
n  n 

-  k2C  ) 
a* 

P 

P 

C*  =  ik<j>  (C^ 

-  k2cr.) 

P 

p 

Yawing  Motion  (r  = 

4-) 

C  =  ik<J>(C 

-  k2C  ) 

y  y 

y_ 

r 

r 

Cn  =  ikc^(Cn 

-  k2C  .. ) 
n 

r 

r 

CA  =  ik4>(C^ 

r 

-  k2cr. ) 

r 

k  d>  C  +  d>  C 

y*  y, 

p  <t 


k  d>C  +  <+» C 
n*  n , 

P  <t> 


k  *CA.  +  <*C* 


k  iLC  +  d/C 

y.  y 


k  (1/C  +  d/C 

n*  n, 

r  4; 


k  4;CJl.  +  *CX 


* 


(37) 

(38) 

r  39) 


(40) 

(41) 


(42) 


2.8  Dynamic  Stall  Effect 

Dynamic  stall  effect  is  important  in  predicting  dynamic 
stability  parameters  for  configurations  with  low-swept  wing 
planforms.  Since  theoretical  calculations  are  very  time-consuming 
and  may  not  be  accurate,  it  is  proposed  to  calculate  three- 
dimensional  results  by  using  2-D  test  data.  The  concept  is  quite 
similar  to  that  in  steady  flow  to  account  for  the  effect  of  boundary 
layer  separation  (ref.  1).  That  is,  iterative  calculation  is  first 
made  for  steady  separation  effect  to  obtain  local  effective  angles 
of  attack  (ae).  Local  effective  oscillation  amplitudes  are 
calculated  by  using  a  three-dimensional  unsteady  attached-flow 
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method  (ref.  4).  Two-dimensional  test  data  are  interpolated  based 
on  these  local  effective  angles  of  attack  and  oscillation 
amplitudes.  To  interpolate,  the  2-D  test  data  must  be  expressed  in 
terms  of  nonlinear  functions  of  a  ^  and  which  depends  on 
oscillation  amplitudes.  This  is  done  by  Fourier  analysis.  The 
resulting  mean  value  is  a  function  of  mean  angles  of  attack  (am) 
which  are  interpolated  for  the  appropriate  sets  of  data  to  be  used 
at  the  calculated  ag.  One  set  of  calculations  will  be  illustrated 
later  to  expose  this  idea. 


3.  NUMERICAL  RESULTS  AND  DISCUSSION 


In  the  following,  calculated  results  for  delta  wings  with  sweep 
angles  of  60°  and  80°  will  be  presented.  In  addition,  calculated 
results  for  the  F-106B  configuration  will  be  compared  with  data. 

All  calculations  are  made  in  eight  iterations  to  adjust  free  vortex 
positions.  The  vortex  breakdown  effect  is  applied  only  at  the  last 
iteration.  Finally,  calculated  results  for  a  rectangular  wing  with 
dynamic  stall  effect  will  be  illustrated. 

Delta  Wing  with  a  Sweep  Angle  of  60° 

Lateral  coefficients  for  this  configuration  are  presented  in 

Fig.  1.  The  aerodynamic  characteristics  are  significantly 

influenced  by  free  vortex  position,  the  motion  lag,  and  the  bursting 

locations.  Typically,  C  is  positive  without  vortex  breakdown  and 

becomes  negative  only  with  vortex  bursting.  As  it  can  be  seen,  due 

to  early  vortex  breakdown  for  this  configuration,  the  resulting 

rolling  moment  coefficients  due  to  p  become  largely  negative  and 

cause  the  roll  damping  derivatives  to  be  more  negative.  Also  as  a 

is  increased,  the  effect  of  j3  derivatives  due  to  lag  in  motion  is 

more  appreciable.  As  the  reduced  frequency  is  increased,  the 

magnitude  of  C  will  diminish  but  still  is  negative.  The  test  data 

0 

indicate  that  vortex  bursting  occurred  at  a  lower  a  than  the 
calculated  value  which  was  based  on  a  steady-flow  model.  It  should 
also  be  noted  that  the  vortex  bursting  point  is  known  to  be  affected 
by  dynamic  motion.  Since  a  systematic  set  of  data  is  not  available, 
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the  dynamic  effect  on  vortex  bursting  cannot  be  assessed  at  the 
present  time. 

The  discrepancy  between  the  calculated  results  and  experiments 
is  mainly  due  to  prediction  of  vortex  breakdown  location. 

Delta  Wing  with  a  Sweep  Angle  of  80° 

Lateral-directional  derivatives  for  this  configuration  are 

determined  at  a  reduced  frequency  of  0.06  and  an  amplitude  of  5°. 

The  calculated  results  show  large  positive  rolling  moment 

coefficients  due  to  (3  (see  Fig.  2).  For  this  configuration,  the 

vortex  breakdown  occurs  at  about  a  =  38°.  Therefore,  this  positive 

rolling  moment  due  to  {3  causes  to  be  more  positive,  which  may 

P 

result  in  wing  rock.  The  correlation  between  the  experimental  data 
and  numerical  results  is  good.  At  high  angles  of  attack,  the  effect 
of  (3  derivatives  is  more  significant  and  is  therefore  important  to 
dynamic  stability  of  airplanes,  in  particular  those  with  highly 
swept  wings.  The  present  numerical  results  show  that  the 
acceleration  derivatives  can  be  predicted  within  a  reasonable 
accuracy  for  highly  swept  wings.  Also,  in  a  forced-oscillation 
test,  it  is  not  possible  to  decouple  the  effect  of  (3  from  the  P 
motion.  However,  by  the  present  method,  it  is  possible  to  decouple 
the  effect  of  each  mode  separately. 
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F-106B  Configuration 

This  configuration  has  a  conically  cambered  wing  of  60-deg 

sweep  and  a  vertical  tail.  At  a  reduced  frequency  of  0.2,  the 

factor  for  the  remaining  vortex  lift  after  breakdown  is  not  applied, 

since  not  enough  data  on  dynamic  vortex  breakdown  are  available  for 

guidance.  The  calculated  results  are  presented  in  Fig.  3(a)  for 

roll  stability  derivatives  and  in  Fig.  3(b)  for  yaw  derivatives.  It 

is  seen  from  Fig.  3  that  predicted  roll  and  yawing  moment 

derivatives  are  significantly  improved  by  including  {^-derivatives . 

Note  that  C  at  a  =  35  deg  is  quite  negative  and  is  contributed  by 
y.  -  - 

the  vertical  tail.  Since  the  dynamic  fuselage  effect  is  not 
available  in  the  current  version  of  the  code,  it  is  not  known 
whether  it  will  improve  the  predicted  results  at  a  =  35  deg. 

Rectangular  Wing  of  Aspect  Ratio  4  with  Dynamic  Stall  Effect 

Some  dynamic  stall  data  on  several  airfoil  sections  are 

available  in  reference  11.  However,  only  data  for  the  NACA  0012 

airfoil  are  more  complete  and  therefore  will  be  used  in  the 

following.  This  set  of  data  was  obtained  by  oscillating  the  airfoil 

about  the  quarter-chord  point  at  mean  angles  of  attack  (a  )  of  3,  4, 

m 

5,  and  10  degs  and  an  amplitude  of  10  deg  at  a  reduced  frequency 

(0)0/2^)  =  0.1  and  a  Mach  number  of  0.3.  The  sectional  c^  -  a 

curves  are  presented  in  Figure  4.  These  sectional  data  are  then 

analyzed  through  Fourier  analysis  to  obtain  c^  as  functions  of 

and  a  .  For  example,  at  a  =5  deg,  c,  can  be  written  as 
1  in 
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cz  =  0.4825  +  (5.114  +  1.184^  +  5.198a^  +  79.694aJ 
+  42.474aJ  -  0.005iJ  -  0.2298a1aJ  -  0.245a^a^)a1 
+  (0.007  +  0.0979^  -  0.051a^  -  2.171aJ  -  7.149aJ 

+  0.006(1^  +  0.0412aJaJ  -  0.001^  -  0.0328aja1)a1  (43) 

where  the  constant  term  on  the  right-hand  side  is  the  mean  value 

of  (i.e.  ).  In  the  calculation  for  a  wing,  this  c^  is 

m  m 

linearly  interpolated  with  the  effective  angle  of  attack.  In  Figure 

5,  the  calculated  aerodynamic  models  at  ctm  equal  to  5  and  10  deg  are 

compared  with  data.  The  agreement  is  reasonably  good. 

Assume  that  a  rectangular  wing  of  aspect  ratio  4  and  with  the 

NACA  0012  airfoil  section  is  oscillated  in  pitch  about  the  quarter 

chord  point  with  the  same  flow  conditions  as  those  used  in  the 

airfoil  testing.  The  resulting  total  -  a  curve  is  plotted  in 

Figure  6  and  calculated  sectional  c^'s  at  three  spanwise  stations 

are  illustrated  in  Figure  7.  Since  there  are  no  experimental  data 

for  the  wing  available  for  comparison,  no  assessment  of  accuracy  can 

be  offered.  However,  the  results  appear  to  be  plausible.  Note  that 

calculated  c^  and  do  not  pass  through  the  coordinate  origin. 

This  discrepancy  is  caused  by  the  assumption  of  linear  interpolation 

of  c^  with  respect  to  the  effective  angle  of  attack.  Except  having 
m 

this  discrepancy,  Eq.  (43)  and  other  similar  equations  not  listed 
here  appear  to  be  capable  of  modeling  the  dynamic  stall  effect  at 
different  amplitudes  and  mean  angles  of  attack. 
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In  applications  to  fighter  configurations,  dynamic  stall  data 
of  typical  thin  airfoils  are  needed.  These  data  can  then  be  used  to 


predict  dynamic  effect  on  longitudinal  and  lateral-directional 
aerodynamic  characteristics. 
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4.  CONCLUSIONS 


An  unsteady,  compressible,  aerodynamic  method  was  developed  to 
calculate  dynamic  stability  parameters  at  high  angles  of  attack.  In 
the  method,  the  edge-separated  vortex  flow  was  represented  by  free- 
voitex  filaments  in  frequency  domain.  The  main  emphasis  of  the 
method  was  to  calculate  dynamic  lateral-directional  stability 
parameters.  Calculated  results  of  C  for  a  delta  wing  of  80  degs 

P 

were  shown  to  produce  total  roll  damping  derivatives  which  agreed 

well  with  data.  For  the  60-deg  delta  wing,  the  calculated  trend 

of  C  variation  with  a  is  correct,  except  that  its  magnitude  is 
f3 

lower  at  low  a  than  the  data  showed.  Calculated  results  for  an  F— 
106B  configuration  with  the  unsteady  method  show  significant 
improvement  in  predicting  dynamic  roll  and  yaw  stability 
derivatives. 

Calculation  of  dynamic  stall  effect  for  a  rectangular  wing  of 
aspect  ratio  4  was  also  demonstrated.  The  calculation  was  based  on 
interpolation  of  two-dimensional  test  data  at  the  calculated 
sectional  effective  angles  of  attack.  No  3-D  data  were  available 
for  comparison.  However,  the  calculated  results  appeared  plausible. 
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5.  RECOMMENDATIONS 


This  investigation  represents  an  initial  attempt  to  predict 
dynamic  effect  on  aerodynamic  characteristics  at  high  angles  of 
attack.  To  improve  the  prediction  method,  the  following 
experimental  and  theoretical  work  would  be  needed. 

(1)  Obtain  experimentally  a  systematic  set  of  dynamic  stall 
data  in  c^,  cd,  cm  on  a  thin  airfoil  typical  of  current 
fighter  wings.  This  set  of  data  can  be  used  in  the 
present  method  to  calculate  dynamic  stability  derivatives 
on  wings  with  low  sweep  angle. 

(2)  Obtain  experimentally  a  systematic  set  of  dynamic  vortex 
breakdown  data  on  delta  wings  of  different  sweep  angles. 
The  data  should  include  the  dynamic  effect  on  bursting 
position  and  total  lift  variation.  This  set  of  data  will 
be  analyzed  and  incorporated  into  the  present  method  to 
calculate  dynamic  stability  derivatives  on  wings  with  high 
sweep  angle. 

(3)  Develop  an  unsteady  mathematical  model  to  represent  the 
fuselage  effect. 
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10  15  20 

a,  deg. 


J 5  40 


Figure  1  Rolling  Moment  Coefficient  due  to  Side  Acceleration 
for  a  60-deg  Delta  Wing  at  M  =  0.13,  ’k  =  0.066  and 
Amplitude  =  5  deg. 
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O  Data  (ref.  13) 
Theory 

-  Unstead  Flow 

- - —  Steady  Flow 


0  Data  (ref.  14) 


Figure  3  Dynamic  Stability  Derivatives  for  an  F-106B 
Configuration  at  M  =  0.1  and  k  =  0.2.  Moment 
Center  at  G.275  c. 
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Figure  5  Aerodynamic  Models  for  an  NACA  0012  Airfoil 
in  Dynamic  Stall  Obtained  from  Figure  4. 


0.0  2.0  4.0  S.O  3.0  1  0.0  12 .0  14.0  !  6.0  13.0  2 0.0  22.0 

a,  deg. 

Figure  6  Calculated  Lift  Response  on  a  Rectangular  Wing  of 

Aspect  Ratio  4  with  NACA  0012  Airfoil  in  a  Pitching 
Oscillation  about  l/4c  at  M  -  0.3  anf  k  *  0.1. 
a  =  10°  +  10°  sinut. 
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APPENDIX  A: 


DERIVATION  OF  INDUCED  VELOCITY  DUE  TO  AN 
OSCILLATING  LINE  VORTEX 


The  expression  for  induced  velocities  due  to  an  oscillating 
horseshoe  vortex  are  available.  To  obtain  the  corresponding 
expressions  for  an  oscillating  line  vortex  in  compressible  flow,  the 
effect  of  trailing  legs  in  the  horseshoe  vortex  must  be  removed.  To 
achieve  this  objective,  the  general  equation  for  oscillating 
doublets  is  needed  (ref.  4): 
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By  integration  by  parts,  equation  (A. 2)  becomes 
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x0=x-^=x-x1-  x(x2  -  Xj) 

y0=y_Ti=!y“yi_  ^  Cy2  -  yp 

The  straight  line  L  is  now  defined  by  (0,  1)  in  x  (see  Fig.  A. 1). 

If  F  in  the  first  term  of  equation  (A. 4)  is  replaced  by  the 

following  form,  equation  (A. 4)  is  reduced  to  oscillating  line  vortex: 

x  2 (x .  -  x) 

F  *  /  ( — ~)z  dr)  +  Arctan[ - ---■  —  . . . . .  . ] 

L  Sr:  °  ,  w,  ,2  .  2 ,  s 2  .  02  2. 
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-  Arctan  {- 
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' . . . - .  . - . - - - 1 

(y9  -  yV(x9  -  x)2  +  £2(y7  -  y)2  +  p2z2 


(A. 6) 


is  differentiated  with  respect  to  z  with  equation  (A. 6)  for  F  and 
letting  the  frequency  go  to  zero,  the  exact  steady  expression  for  a 
line  vortex  element  can  be  obtained  as  will  be  shown  later.  Now 
consider  the  following  expression  in  $ i.e., 
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Let  x^  or  X2  in  equation  (A. 7)  go  to  infinity;  hence: 
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The  lower  limit  of  integration  goes  to  positive  infinity  as  x^  or  X2 
goes  to  infinity;  therefore, 


Aim  1=0 

X.+® 

X2>a> 


Therefore,  it  can  be  concluded  that  ^  is  due  t0  vortex  shedding 
only  and  comes  only  from  bound  vortex  of  oscillating  horseshoe 
vortex.  Thus,  the  equation  for  an  oscillating  line  vortex  becomes: 


£--♦!+  *2  (A>10) 


where  and  ^  are  given  by  equations  (A. 6)  and  (A. 3). 

To  derive  expressions  for  induced  velocities  due  to  an 
oscillating  line  vortex,  equation  (A. 10)  is  differentiated  with 
respect  to  x,  y,  and  z,  respectively. 
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Expressions  from  <p  ^ : 
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Taking  the  derivative  of  <J>  ^  with  respect  to  z  results  in 
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Q  =■  (x2  -  xL)(y  -  yL)  -  (x  -  xL)(y2  -  y^ 

By  applying  the  following  transformation,  the  line  integrals  in 
equation  (A. 14)  can  be  reduced  to  a  definite  integral. 
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The  second  integral  in  equation  (A. 14)  becomes: 
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The  last  integral  in  equation  (A. 14)  is 
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where : 

a  =*  (x2~x1  )2  +  P2)y2~y1)2 
b  =■  -2[(x-x1)(x2-x1)  +  p2(y-y1)(y2-yI)] 
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To  evaluate  these  integrals,  parabolic  approximation  of 
following  form  is  used: 
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Expressions  from  $2 

Applying  integration  by  part  to  ,  the  result  becomes 
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Differentiating  $2  witl1  respect  to  z  gives 
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To  evaluate  the  integrals  in  the  term,  the  following 

approximation  may  be  used  in  I,  1^,  and  I2  terms  (ref.  15): 
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Therefore,  by  applying  the  above  approximation  (eq.  A. 33),  I,  I 


I2  can  be  integrated  exactly.  Note  that  integrals  in  the  form  of 
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can  be  reduced  into  form  like  I  by  taking  integration  by  parts.  For 
P,  twice  integration  by  part  is  necessary. 

By  the  same  approach,  the  u  and  v  veloctiy  can  be  obtained. 

To  check  the  expressions,  it  should  be  noted  that  the  induced 
downwash  due  to  an  oscillating  line  vortex  can  be  reduced  to  the 
downwash  velocity  due  to  a  steady  line  vortex  if  the  oscillating 
frequency  (co )  is  set  to  zero.  In  this  case,  equation  (A.  14)  becomes 


d£  |  =  dF  | 

dz  '  u)=0  8u  dz  '  L 


(A. 35) 


where : 


x 

F-J  (— ^Zgd-ri+Tan-1  [■ 


z(x^-x) 


L  Rr 


-] 


(y1-yV(x1x)2+p2(y1y)2+p2*2 


-  Tan-1  [- 


z(x2~x) 


-] 


(y2-y)/(x2-x)2+T(y2-y)2+j32z2 


(A. 36) 


and 


dK  L _ 1 _  * 

dz  L  o  n  o  y 

/(x-x2)z+r(  y-y2r+r(z-z2r 

(x2-x1)(x2-x)+p2(y2-y1)(y2~y)+p2(z2-2)(z2-z1)^ 

q"+(z-z2)2[ (x2-x1)2+^(y2~y1)2] 


A- 11 


1 


✓(XTc1)2+p2(y-y1)2+p2(z-Zl)2 


(x2-Xj )(x1-x)+p2(y2-y1)(y1-y)+p2(z1-z)(z2-z) 

{  Q  O  O  ""  O  O  O  ^ 


Q2+(z-z1)2[(x2-x1)2Hf  2(y2-y1)2] 


(A. 37) 


I' 


Equation  (A. 37)  exactly  matches  the  downwash  velocity  due  to  a 
steady  line  vortex. 
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